Comparison of structure and transport properties of concentrated hard and soft 

sphere fluids 
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Using Newtonian and Brownian dynamics simulations, tlie structural and transport properties 
of hard and soft spheres have been studied. The soft spheres were modeled using inverse power 
potentials [V ~ with 1/n the potential softness). Although, at constant density, the pressure, 
diffusion coefficient and viscosity depend on the particle softness up to extremely high values of 
n, we show that scaling the density with the freezing point for every system effectively collapses 
these parameters for n > 18 (including hard spheres), for large densities. At the freezing points, the 
long range structure of all systems is identical, when length is measured in units of the interparticle 
distance, but differences appear at short distances (due to the different shape of the interaction 
potential). This translates into differences at short times in the velocity and stress autocorrelation 
functions, although they concur to give the same value of the corresponding transport coefficient (for 
the same density to freezing ratio) ; the microscopic dynamics also affects the short time behaviour of 
the correlation functions and absolute values of the transport coefficients, but the same scaling with 
the freezing density works for Newtonian or Brownian dynamics. For hard spheres, the short time 
behaviour of the stress autocorrelation function has been studied in detail, confirming quantitatively 
the theoretical forms derived for it. 



I. INTRODUCTION 

Monodisperse spherical particles with short range re- 
pulsions show a simple phase diagram, with only one fluid 
at low density and one crystal phase at high density, and 
a first order transition in between^. The paradigmatic 
case is the system of hard spheres (HS), where the in- 
teraction potential is infinite whenever particles overlap 
and zero elsewhere, and therefore there is no energy scale 
in it. In HS, the crystallization transition was first rec- 
ognized by Alder and Wainwright^, and has the freezing 
and melting points at volume fractions 0/ = 0.494 and 
(j)m — 0.545, respectively^. Additionally, a good approxi- 
mation for the equation of state for HS fluid was given by 
Carnahan and Starling using the virial expansion^, tested 
experimentally using sedimentation of screened charged 
colloids^. For slightly soft spheres (SS), where interpene- 
tration is not yet an issue, the phase diagram changes 
quantitatively, shifting the fluid-solid coexistence, and 
introducing temperature as a new variable (the energy 
scale, e, is set by the potential, and e/fc^T, where fc^T 
is the thermal energy, is coupled to the volume frac- 
tion). Inverse-power potentials {V{r) ~ r~") have been 
widely used to model SS, and the phase diagram stud- 
ied as a function of the softness by simulations^ii and 
theories®'^. Because of the qualitative similarity of the 
phase diagram of SS with that of HS, a mapping of the 
former to the latter has also been tried by means of effec- 
tive diametersiSiiiii^. However, an interesting pathology 
arises in HS which is absent in SS: the elastic modulus 
for large frequencies diverges in HS ^'^'^'^i^^ . 

Experimentally, different colloidal systems have been 
used to model HS, such as silica^S or latex particle o^'^^ , 
or collapsed microgeli^ or core-shelliiS, particles. Because 
the hard sphere interaction is an idealization of quasi- 



rigid spheres, these systems are usually referred to as 
" nearly hard spheres" , and the question of their softness 
arises naturally—. In most cases, additionally, a (gen- 
erally thin) polymer layer is adsorbed onto the particles 
to provide steric stability, what increases the softness of 
the particles^*'. Even more, experiments with really soft 
particles, such as swollen microgels2i or star polymers^-, 
can show new phases due to interpenetration, but also 
phase diagrams qualitatively similar to the HS one for 
the harder-particle limit. However, the short range in- 
teraction is not directly accessible experimentally^^, and 
exact knowledge of the true softness of the particle is 
missing. It is therefore important to know what proper- 
ties depend on the softness of the particle, and if there is 
a simple scaling for different systems with different soft- 
ness. 

Recent simulations of inverse power potentials have 
shown that the transport coefficients depend strongly on 
the particle softness, and the HS limit is reached for large 
values of n (n > 72)2^. It has been shown that the time- 
correlation functions (and transport coefficients) can be 
approximately scaled using the exponent v^^^, but nei- 
ther the scaling is perfect nor the exponent is known 
experimentally. On the other hand, the transport coeffi- 
cients behave quasi-linearly with the inverse packing frac- 
tion for different softness, what can indicate that the key 
parameter for mapping the systems is the free volume, 
but a detailed analysis of the parameters shows inconsis- 
tencies (for instance, the free volume vanishes at different 
density for every transport coefficient It is desirable 
to find a simple mapping from soft particles to HS, if it 
exists, that can be applied in experimental systems. 

The freezing point, on the other hand, is well defined 
fundamentally, and can be identified with several well- 
known, albeit approximate, criteria-'^, independently of 
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the interaction potential, such as the Hansen- Verlet cri- 
terion, which uses the height of the neighbour peak in 
the structure factor—, and a dynamical one, proposed 
by Lowen and co-workers, based on the decrease of the 
long time diffusion coefficients^. In computer simula- 
tions, the freezing point can be determined accurately 
using Gibbs-Duhem integration from HS-. Therefore, in 
this work, we propose to use the density at freezing to 
scale the results. Using computer simulations, we show 
that the structural and dynamical quantities of interest 
of inverse-power potentials can be rationalized when plot- 
ted against the scaled density, p/pfroezing, for large den- 
sity at identical temperature. The long range structure 
is almost identical for systems with equal /o/yOfreezing, but 
differences appear at short distances due to the differ- 
ent interaction potentials. The pressure, which shows a 
rather complicated trend when studied at constant den- 
sity, increases monotonously with the softness at constant 
p/pfrcozing- The diffusion coefhcient and viscosity can be 
collapsed for n > 18, including HS, and the time correla- 
tion functions only show differences at short times, giving 
different shear moduli at large frequencies. Our conclu- 
sion is therefore that the relative density to the freezing 
transition is the key parameter governing the structure 
and dynamics of the system, for large enough n and den- 
sity. This poses a simple criterion that can be useful for 
experiments. 



which we apply standard periodic boundary conditions. 
Simulations are started from a simple cubic lattice. Ve- 
locities are then drawn at random from a Gaussian dis- 
tribution to let the system evolve towards equilibrium. 
From such an equilibrated configuration all further sim- 
ulations are run. 

The standard deviation of the velocity determines the 
temperature of hard spheres: (f^) = -^ksT. Inherently 
Newtonian dynamics does not show any dependence on 
temperature. Changing the temperature leads to a dif- 
ferent timescale, but does not modify the dynamics qual- 
itatively. It is therefor the density alone that determines 
the physical behaviour of the hard spheres. 

In all our simulations the diameter a determines the 
length scale and the particle mass m the unit of mass, 
while fcsT may vary. With these conventions the time is 

given in units of r = </ \J T^T^ wherever it appears. 

We have introduced Brownian motion by incorporat- 
ing a random thermostat. This feature introduces a new 
timescale tb on which the random kicks of the solvent 
are mimicked by drawing new velocities for each parti- 
cle from a Gaussian distribution^*^. This timescale must 
be chosen such that it effectively generates diffusion be- 
tween collisions. Given densities up to the freezing point 
we find that tb = 0.017 is sufficiently small (smaller val- 
ues of Tb cause too slow microscopic dynamics, increasing 
the simulation time^i). 



II. SIMULATION DETAILS 

We make simulations of different monodisperse sys- 
tems: hard spheres and soft spheres with different "soft- 
ness". Whereas the interaction potential is continuous 
for soft spheres, it is not for hard spheres, what makes an 
important difference in the computational method. The 
microscopic dynamics is Newtonian, but we have also 
performed some Brownian dynamics or damped Newto- 
nian dynamics to make more direct contact with experi- 
mental colloids. 



A. Hard spheres 

Hard spheres only interact upon contact, when they 
collide elastically due to their excluded volume. Given 
such a discontinuous potential we resort to an event- 
driven algorithm for the hard-sphere simulations. The 
method relies on the fact that in between collisions par- 
ticles move undisturbed and in a straight line. It is there- 
fore advantageous to only process the collisions when 
they occur. This is the basic feature of an event-driven 
algorithm: it docs not evolve with a fixed time step, but 
instead time propagation is determined by the events in 
the simulation — the hard sphere collisions in the case 
of Newtonian dynamics. In this paper we have modeled 
the interactions of the hard spheres by elastic collisions. 

The simulation consists of 1000 particles in a cube to 



B. Soft spheres 

For soft spheres we use a continuous repulsive poten- 
tial: 

Usc{r) = e (-) (1) 

where e sets the energy scale and s = 1/n is the softness 
of the potential. In this work, we change n from n = 10 
to n = 36, what provides a wide enough range to show 
clear effects from the softness of the potential, and the 
crossover to effective hard spheres. Figure [1] presents the 
interaction potentials for n = 36, 24, 18, 12 and 10 and 
the hard core repulsion. Note that as n decreases (the 
softness increases), the potential has a longer range, and 
for n = 10, the interaction energy is negligible only for 
distances above r « 2cr. 

Since the interaction potential is soft, the particle vol- 
ume cannot be unambiguously defined, and thus the vol- 
ume fraction is not a convenient control parameter for 
these systems, although it is commonly used in experi- 
ments. Instead, states in potentials with different soft- 
ness can be compared using the effective density defined 
by the parameter F — pa^{e/kBTY^", where p is the 
number density. In this work, the temperature is fixed to 
fc^T = e = 1, and variations in F are caused by changes 
in the number density exclusively. 
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Figure 1: Pair interaction potentials used in this work, with 
decreasing n (increasing softness) from below: n = 36, 24, 
18, 12 and 10, and hard sphere. Note the long range of the 
repulsion for the lowest n. The inset shows the freezing and 
melting lines in terms of F— . The stars indicate state points 
checked in simulations. 



Conventional molecular dynamics are performed inte- 
grating the Newton's equations of motion with veloc- 
ity rescaling to maintain constant temperature. The 
time step f or th e integration is 5t = 0.002, in standard 
units of a\pmjt. All potentials have been truncated for 
V < 10-3. 

Damped Newtonian dynamics have been performed 
also, introducing friction and fluctuation (Brownian) 
forces. For particle j, the equation of motion is now 
Langevin equation: 
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0.970 


18 


0.0555556 


1.016 


12 


0.0833333 


1.167 


10 


0.10 


1.300 



Table I: Tabulated data of the equilibrium phase diagram 
from Agrawal and Kofke^. 




Figure 2: Pressure of the different systems studied here, 
as labeled, and Carnahan-Starling equation (magenta line). 
The vertical dashed lines indicate the freezing densities from 
Agrawal and Kofke^. 



i 

where Fij is the interaction force between particles i 
and j, 7 is the friction coefficient with the solvent, and 
fj is the Brownian force. The fluctuating and viscous 
terms are related by the fluctuation dissipation theo- 
rem: {fj{t')fi{t)) = 6kBTjS^j6{t - t'). In this work, 
7 — 2Q^/em/a which gives a strong damping (the mo- 
mentum relaxation time scale, m/7 0.05, is slightly 
smaller than the typical collision time between particles). 
The equations of motions are integrated using the algo- 
rithm developed by Heun^^ , with a time step 5t = 0.0005. 

III. RESULTS AND DISCUSSION 

Agrawal and Kofke used Gibbs-Duhem integration to 
compute the solid-fluid phase boundary for SS modeled 
with inverse-power potentials, starting from the FCC- 
fluid coexistence in HS^. Their results for the freezing 
transition are presented in Table I (and in the inset to 
Fig. [T]) for different values of n and the HS system. Note 
that freezing density does not evolve monotonously with 
the potential softness (at constant temperature ksT — e. 



r = pcr^): upon decreasing n from infinite (increasing 
softness starting from HS), Fficczing hi'st decreases, and 
then increases for n < 36. 

Figure [2] presents the pressure in the fluid phase up 
to the freezing point for different values of n. Whereas 
at low density the pressure is almost independent of the 
softness, important differences are noticed at high den- 
sity, close to freezing. In particular, every system has a 
different pressure at flxed density, or at its own freezing 
point. Also, at low density (not shown in the figure) the 
pressure increases weakly with the softness, whereas at 
high ones it decreases. 

Contrary to other works, we will not try to rationalise 
the complex picture presented in Fig. [2]with the so-called 
effective diameter-^^^, mapping the pressure to the hard- 
sphere one, or other theoretical approximations^, but we 
will show that a convincing scaling of the data can be 
performed using the freezing density, showing that the 
hard sphere limit can be reached, within this scaling, 
for n as low as n « 24. In the upper panel of Fig. [3] 
we present the pressure as a function of the scaled F. 
A monotonous trend of the pressure with the potential 
softness is now evident at all shown densities. 

In the lower panel the pressure has been scaled with 
the value at freezing, showing almost perfect collapse of 
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freezing 

Figure 3: Upper panel: Pressure as a function of the 
scaled density for different potential softness, as labeled, and 
Carnahan-Starling equation (magenta line). Lower panel: 
Scaled pressure vs. scaled density for the same values of n. 
The inset shows the scaling factor for the pressure (red circles) 
and the ratio of second virial coefficients to the hard sphere 
system (blue asterisks) for different n. 



the data for n > 18. The inset presents the scaling factor. 
An, (ratio of the pressure at freezing for SS to pressure for 
HS) and the ratio of the second virial coefficient, B2, to 
the hard sphere value, as a function of the potential 

softness. These two quantities disagree, showing that a 
mapping based on the second virial coefficient would not 
find such a good scaling, concluding that soft spheres 
behave as hard ones only for very large values of n. The 
a priori surprising finding emerges, that a scaling of soft 
spheres onto hard ones works better at high densities, 
where the closeness to freezing provides a common scaling 
variable, which collapses the data much better than the 
second virial coefficient would. 

The scaling observed in Fig. O in contrast with the 
complex situation observed in Fig. [21 indicates that the 
relevant length scale for thermodynamic properties is not 
a, but the average particle separation, r~^/^ = p~^l^ ja. 
This can be observed in the structure of the system, 
which can be scaled for different states and systems when 
p~'^l^ is used as the unit of distance. In Fig. [4] the pair 
distribution function is presented for the freezing density 
at every softness. The inset shows the raw data with 
distance measured in units of cr, where the maxima are 




Figure 4: Pair distribution function of the soft and hard 
spheres at freezing (see table 1), for different softness: from 
top to bottom in the nearest neighbour peak HS, n = 36, 
n = 24, n = 18, n = 12 and n = 10. The inset shows the raw 
distribution as a function of the average distance measured 
in units of a, and the main panel in units of p~^^^ (the hard 
sphere system has been left unchanged). 




Figure 5: Structure factor of the soft and hard spheres at 
freezing (see table I), with the same colour code as Fig. 3] 
The wave vector has been scaled with the inverse average 
distance p^^'^ (the hard sphere S(g) has been left unchanged). 
The inset shows the neighbour peak in detail. 



shifted to lower distances for higher n, and the main panel 
presents g(r) with the distance scaled with the mean in- 
terparticle distance, where scaling has been achieved for 
all peaks (in order to keep the first neighbour peak in 
HS at r = cr, distances have been scaled not simply by 
but by r7^[/^ = where F^^ denotes 

the hard sphere value of F at freezing). Differences in the 
local structure are still observed in the first neighbour 
peak, caused by the differences in the potential softness, 
and these are responsible for the differences observed in 
the pressure (Fig. [3]- upper panel). 

The relevance of the interparticle distance in contrast 
to a for structural properties, can also be shown using 
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the structure factor, which is more accessible experimen- 
tally than the pair distribution function. Fig. [5] presents 
the structure, S{q), as a function of the wave vector 
scaled with the inverse separation distance p^^^ at freez- 
ing for different potential softness. The inset shows the 
neighbour peak in more detail, where the differences in 
the local structure should be noticed. As expected from 
the Hansen- Verlet criterion^^, the height of the peak is 
around 2.85 in all cases, but no difference in the shape of 
the curves or position of the maximum can be observed 
beyond the noise level. This shows that all systems in 
the range of softness studied have similar liquid struc- 
ture (except locally) at their freezing points, and that 
the scaling factor is the mean interparticle distance, thus 
establishing the density that must be used to match the 
thermodynamic quantities in systems with different soft- 
nesses. 

We move now to the study of dynamic quantities, par- 
ticularly the transport coefhcients (diffusion coefhcient 
and the shear viscosity), and the autocorrelation func- 
tions involved in their calculations. The diffusion co- 
efficient, D, can be obtained via the Green-Kubo rela- 
tion from the velocity autocorrelation function, or from 
the long time slope of the mean squared displacement^ii, 

D = r C,{t)dt = hm ^^^^ (3) 

where C^{t) = 1/3 {vi{t) ■ ?7,(0)) while 5r = r{t) - f(0), 
and the brackets indicate ensemble and time-origin aver- 
aging. The velocity autocorrelation function is presented 
in Fig. [6] (upper panel), for SS with Newtonian dynam- 
ics (ND) and damped Newtonian dynamics (dND), at 
the freezing point for every system. Cv{t) decays faster 
for dND due to the random forces and friction with the 
solvent, resulting in lower diffusion coefficient than ND. 
(The presence of backwards motion indicates that the 
Brownian forces are not strong enough to provoke the 
decorrelation of velocities before collisions occur). Be- 
cause the local environment of the particles is different 
for the different systems, the Cv{t) functions do not co- 
incide in the decay and at the backscattering minimum, 
which is deeper for softer potentials. These differences 
however concur to give equal diffusion coefficients for dif- 
ferent particle softness, as shown below in Fig. [T] 

The mean squared displacement (MSD) is plotted in 
the lower panel of the figure for SS with ND and dND, 
and for HS with ND and Brownian dynamics (BD). Due 
to the Brownian forces (in dND) or kicks (in BD), the 
MSD grows slower for dND or BD than for ND. (The 
short time dynamics for HS with BD is controlled by the 
arbitrary time interval between Brownian kicks, tb, and 
the MSD in this case has been shifted in time to give the 
same long time self diffusion coefficient). Note that the 
MSD with the same microscopic dynamics collapse inde- 
pendently of n (except for HS with BD), even though the 
density is rather different. The velocity autocorrelation 
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Figure 6: Velocity autocorrelation function (upper panel) and 
mean squared displacement (lower panel) for SS with n — 36 
(red), n = 24 (green), n = 18 (blue) and n = 10 (brown) and 
HS (only the MSD) (thin black line) at the freezing density 
(see table I) with ND (upper set of curves in both panels) and 
BD or dND (lower set of curves). The Cv{t) for dND have 
been displaced vertically for clarity. BD has been rescaled on 
the time axis by a factor td = 0.1132 such that the long time 
diffusion agrees. 



functions did not collapse because this function focuses 
on the microscopic dynamics where the local structure is 
sampled, although the final diffusion coefficient is iden- 
tical for all potentials. The scaling in MSD, and the 
concomitant agreement of diffusion coefficients, indicates 
again that the relevant quantity for this transport coef- 
ficient is the distance to freezing, instead of the density. 

The diffusion coefficient from the long time limit of the 
MSD with ND for the different potentials is presented in 
Fig. [71 as a function of the density (upper panel), and 
as a function of r/Ffroozing (lower panel). Whereas the 
raw data shows differences due to the softness of the po- 
tential up to n = 36 (upper panel), in agreement with 
previous results^^, the value of the diffusion coefficient 
at the freezing point is very similar in all cases, as ex- 
pected from the dynamic criterion for freezing proposed 
by Lowen et alM. (note, however, that a short-time dif- 
fusion coefficient cannot be defined for ND). When the 
diffusion coefficient is plotted against F/Fficczing, a con- 
vincing collapse of the data is found for n > 24, with 
small deviations for n = 18 (lower panel). The same 
scaling that was used to rationalise the structural and 
thermodynamic properties is successful in the collapse of 
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freezing 



Figure 7: Difltusion coefficient for increasing density approacli- 
ing the freezing transition for soft spheres with different soft- 
ness and hard spheres, as a function of the density, p (upper 
panel) and the scaled density, r/Ffrcozing (lower panel). From 
top to bottom in the upper panel (bottom to top in the lower 
one): n — 10, n — 12, n — 18, n — 24, n = 36 and hard 
spheres. 



the difFusion coelRcient. Similar collapse of the diffusion 
coefficient can be found using the BD or dND as short 
time dynamics. 

Finally, we study the viscosity of the system, which can 
be calculated from the stress autocorrelation function, 
Co-a-{t). The stress tensor is calculated as: 



N 

i=l 



N 

i<j 



(4) 



where Via is the a-th component of the velocity of particle 
i, and Fijp is the /3-th component of the interaction force 
between particles i and j. 

For HS, however, the interaction potential is not dif- 
ferentiable, and the force cannot be derived from it, al- 
though there are collisions between particles. The in- 
teraction force can be found using the exchange of mo- 
mentum between colliding particlesS-i^. To this end 
one needs to average the change of momentum Fijp = 
mAviji36{t — tij) during a brief time interval Tq-, where 
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Figure 8: Stress correlation function at the freezing density 
for different potential softness with Newtonian dynamics (up- 
per panel) and Brownian dynamics (lower panel); from top 
to bottom n = 36, n = 24, n = 18, n = 10 and hard spheres 
(black). Dashed black line in the upper panel refers to sim- 
ulations with To- = 0.0017 to probe shorter correlation times. 
The inset shows the convergence to l/\/i (magenta dotted) for 
sufficiently short time. To probe these timescales a Brownian 
time step tb = 0.0017 was used (dashed line, inset). 



the (5-function refers to the instantaneous change of mo- 
mentum, which leads to 



N 

1=1 



1 

— ^ mrijaAvijp, 



(5) 



where the second summation now refers to all collisions 
within the time interval t^. Unless otherwise stated — 
0.017. 

The stress correlation function, Caa{t) = 
Q<^(o-Q/3(i)o-a/3(0)), IS presented m Fig. [5] for 
the hard and soft spheres with different softness at their 
freezing points, using ND and BD or dND (the brackets 
indicate average over time origin). In both microscopic 
dynamics, differences are noticed in the correlation 
function at short time, but they collapse within the noise 
level at longer times (compared to the microscopic time 
scale) for n > 24. The curve for HS with BD has been 
scaled in time by the same factor as the MSD (see Fig. 
ini lower panel), and the collapse with SS is observed for 
T > 10"^, as in the MSD case. 

The short time divergence of the HS correlation 
function can be observed in the inset of the bottom 
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Figure 9: Calculation of the viscosity from the squared inte- 
grated stress AA{t) at different densities, indicated by color. 
Solid line with filled symbols show IfS, dashed lines indicate 
soft spheres for n = 36. Thin dashed lines mark the theoret- 
ical result of the 5-contribution in the stress autocorrelation 
function with the weight given in eq. ([8]). 
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C..(i->O) = ^0Va)./l^./f, (6) 

OTa V TT V I 

2 

where Ta = -fjj- with I?o being the short time diffusion 
coefficient. (These short times were only accessible by 
reducing the Brownian time scale tb by a factor of 10.) 

For the case of ND, on the other hand, it has been 
shown22, that the stress correlation function of HS has 
two contributions: a (5-function with weight 775 at the ori- 
gin and a well behaved bounded curve for the remaining 
region. We sample the latter part only. 

Similar to the diffusion coefficient, the viscosity can 
be calculated via a Green-Kubo relation with the stress 
correlation function, or by the long-time limit of the in- 
tegrated stress: 

/■°° PiV 1 

V = dtC{t) = i— lim -(AA(t)2), (7) 

Jq 2 t^oo t 

where AA(t) \ /o Eq</3 and the brackets 

indicate averaging over time origin for this integral, and 
V is the volume of the system. In HS, the (5- function 
in Ccrcr{t) mentioned above gives a contribution to the 
viscosity at i = equal to 775 with^: 

192 / Ftt \ 

where P is the pressure and tjq is the viscosity from ki- 
netic theory ryo — ^\/ "'tf ■ Numerically, the in- 
tegration of the correlation function is noisier and there- 
fore, the viscosity is usually calculated with the integral 



Figure 10: Viscosity of hard and soft spheres from the long 
time slope of the integrated stress, with microscopic ND (up- 
per panel). Different potential softness is studied: black 
points, hard spheres; red ones, n — 36; green ones, n = 24; 
blue ones, n = 18; orange ones, n = 10. The lower panel 
shows the Stokes-Einstein relation Drj, with the same colour 
code, and the dashed green thick lines marks the slip value 
l/(27r) (upper line) and stick one l/(37r) (lower line). 



of the stress tensor (second equality in eq. ([7])), presented 
in Fig. [H]- 10000 independent calculations of 77 have been 
performed and averaged. Note the finite short time limit 
for HS, which agrees with the theoretical value 2'qs, in 
stark contrast to the decay to zero for SS. This window 
of qualitative difference between SS and HS moves to 
shorter times upon increasing n. The viscosity, long time 
plateaus in Fig. [51 are finally presented in Fig. [TO] for 
hard and soft spheres with microscopic ND, as a function 
of r/Ftreozing- The crror bars are the standard deviation 
of the averaging. (Note that Fig. [Hi is in logarithmic 
scale, what enlarges the difference between HS and SS at 
small values of the viscosity.) 

Due to the noise of the data, the collapse is not as im- 
pressive as for previous quantities. Deviations are only 
noted beyond the noise level at high density, close to 
Ffreczing, for Ti = 10. Thcsc rcsults are compatible with 
the previous ones in showing that the appropriate scal- 
ing of the thermodynamic quantities and transport co- 
efficients is obtained when F is referred to its value at 
freezing, although definite proof is only obtained for pre- 
vious quantities. We have included it here, nevertheless, 
because of its relevance to experiments in colloid science. 
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The lower panel of the figure shows the product Drj, 
which would be constant according to the Stokes-Einstein 
relation. Surprisingly, as this relation derived from con- 
tinuum fluid mechanics need not apply to the present sit- 
uation of the motion of a particle surrounded by a fluid 
of identical particles, within the noise level, this prod- 
uct is indeed constant, and coincides with the theoretical 
value of the Stokes' law for a sphere with slippery surface 
Drj = (27rCT)-i for n = 36. 

IV. CONCLUSIONS 

We have studied in this work the structural and trans- 
port properties of hard and soft spheres with different 
softness using simulations with different microscopic dy- 
namics. The pressure and the transport coefficients can 
be rescaled onto the HS values at high density when stud- 
ied as a function of the density normalized to the freez- 
ing point for every system for n > 18, although they 
typically show non-monotonous behaviour for increasing 
density (or volume fraction), irrespectively of the mi- 
croscopic dynamics. Concomintantly, distances must be 
measured in units of the interparticle distance, and all 
systems have the same long range structure at their freez- 
ing points. Differences due to the different softness are 
appreciable only in the local environment of the parti- 
cles, and at short times in the velocity and stress auto- 
correlation functions. The latter is particularly interest- 
ing since its value at zero-time depends on the particle 
softness and diverges for hard spheres. The analytical 
prediction at short times for the stress autocorrelation 
function of hard spheres and the initial ^-contribution to 
the viscosity have been verified quantitatively with our 
simulations. As expected, differences in the correlation 
functions are also noticed depending on the microscopic 



dynamics used, although the final values of the transport 
coefficients can be scaled. 

These results show that the key parameter for dense 
systems of hard and soft spheres is not the density, but 
the relative distance to the freezing point (which itself 
shows a non-monotonous behaviour with the softness of 
the potential). This can be of special interest for exper- 
imentalists, where a small softness cannot be ruled out, 
in particular for particles covered with surfactant layers, 
or in cases where the determination of the colloid den- 
sity poses a major problem. Importantly, and perhaps 
at first counter intuitively, our simulations show that the 
mapping of a fluid of soft spheres onto the one of hard 
spheres becomes progressively better at higher densities 
and closer to freezing. This supports the idea, going back 
to van der Waals, that the local packing of a dense fluid 
is dominated by the rapidly varying repulsion which can 
successfully be modeled by the hard sphere onsi^i. In 
these cases, our results indicate that the relevant quan- 
tity to compare different systems is not the bare density, 
but its ratio to the freezing density. In passing, we also 
note that this scaling is much simpler and more effective 
than other attempts to map the structure and dynamics 
of soft spheres onto those of hard spheres, based on ef- 
fective hard-core diameters (as far as the particle is not 
extremely soft, i.e. n > 18). 
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